################################################################################
## VISUALIZATION SCRIPT FOR PRAGLOOK
## create spaghetti and dwell time plots for praglook experiment
##
## jdegen 6/19, based on mcf 6/13
################################################################################

## PRELIMINARIES
# set your path
setwd("~/cogsci/projects/eye_tracking_code/et-ana/")

## PRELIMINARIES
rm(list = ls())
source("useful.R")
source("et_helper.R")

d <- read.csv("processed_data/praglook processed.csv")

## minor odds and ends
d <- subset(d,stimulus != "cross_white") # remove fixation cross
d$stimulus <- to.n(d$stimulus) # convert to numeric

################ PRELIMINARIES #################
## 1. Read in the orders and merge them with the data
order <- read.csv("info/order1.csv")

nrow(d) # first check number of rows
plot(d$stimulus) # now check the stimulus ordering

# now join in the orders
d <- join(d, order) # use join rather than merge because it doesn't sort

plot(d$stimulus) # check that nothing got messed up
nrow(d) # check the number of rows again

## 2. Define the target ROIs (regions of interest)
rois <- list()
rois[[1]] <- c(0,0,840,550) # left
rois[[2]] <- c(840,0,840,550) # right
rois[[3]] <- c(420,550,840,550) # center
names(rois) <- c("L","R","C")
roi.image(rois)

# use check code to make sure that ROIs look right
d$roi <- roi.check(d,rois) 

# see how the distribution of ROIs looks
qplot(roi,data=d)

# set up correctness
d$correct <- d$roi == d$targ.pos

## 3. Align trials to the onset of the critical word
d <- rezero.trials(d)

## 4. subsample the data so that you get smooth curves
##    I like to do this when I don't have much data so that I'm not distracted 
##    by the variation in the data, but then relax the subsampling if I have more data.
subsample.hz <- 10 # 10 hz is decent, eventually we should set to 30 or 60 hz
d$t.crit.binned <- round(d$t.crit*subsample.hz)/subsample.hz # subsample step


################ ANALYSES #################
# every analysis has two parts: an aggregation step and a plotting step
# - aggregation averages over some kind of unit of interest, e.g. trial type
# - and then plotting is making a picture relative to that aggregation

## 1. SPAGHETTI PLOT OF TRIAL TYPE
# 1. Aggregate your data. Which variables do you need to aggregate by? Tip: look at the plot_spaghetti.R script from yesterday. #
ms <- aggregate(correct ~ t.crit.binned + trial.type, d, mean)

# 2. Make a basic spaghetti plot of target (correct) fixations by trial. 
ggplot(ms, aes(x=t.crit.binned, y=correct, colour=trial.type)) +
  geom_line(size=2)       

# 3. Add a horizontal line indicating chance. Add a vertical line indicating onset of target word.
ggplot(ms, aes(x=t.crit.binned, y=correct, colour=trial.type)) +
  geom_line(size=2) +
  geom_hline(yintercept=.33,lty=2) + 
  geom_vline(xintercept=0,lty=3) 

# 4. Improve the appearance of the plot: rename the x and y axes. Set y axis limits to go from 0 to 1. Set x axis limits to go from -2 to 3. Make the axes start at 0.
ggplot(ms, aes(x=t.crit.binned, y=correct, colour=trial.type)) +
  geom_line(size=2) +
  geom_hline(yintercept=.33,lty=2) + 
  geom_vline(xintercept=0,lty=3) + 
  scale_x_continuous(limits=c(-2,3),expand = c(0,0), name = "Time (s)") + 
  scale_y_continuous(limits=c(0,1),expand = c(0,0), name = "Proportion correct looking") 

# 5. Add error bars with bootstrapped 95% CI
# Get the confidence intervals:
mss <- aggregate(correct ~ t.crit.binned + trial.type + subid, d, mean)
ms <- aggregate(correct ~ t.crit.binned + trial.type, mss, mean)
ms$cih <- aggregate(correct ~ t.crit.binned + trial.type, mss, ci.high)$correct
ms$cil <- aggregate(correct ~ t.crit.binned + trial.type, mss, ci.low)$correct
ms$YMin = ms$correct - ms$cil
ms$YMax = ms$correct + ms$cih
  
# Use either geom_errorbar() or geom_pointrange() to plot error bars
ggplot(ms, aes(x=t.crit.binned, y=correct, colour=trial.type)) +
  geom_line(size=2) +
  geom_hline(yintercept=.33,lty=2) + 
  geom_vline(xintercept=0,lty=3) + 
  geom_errorbar(aes(ymin=YMin, ymax=YMax),position=position_dodge(.05)) +
  scale_x_continuous(limits=c(-2,3),expand = c(0,0), name = "Time (s)") + 
  scale_y_continuous(limits=c(0,1),expand = c(0,0), name = "Proportion correct looking")


## 2. BY ITEM ANALYSIS
# This won't look good until we have a lot of data because we are dividing our 
# data in 6 parts...
# 1. Aggregate your data by time, trial type, and item.
ms <- aggregate(correct ~ t.crit.binned + trial.type + item, d, mean)

# 2. Plot your data as before (without error bars to reduce the clutter), but with different facets for different items.
ggplot(ms, aes(x=t.crit.binned, y=correct, colour=trial.type)) +
  geom_line() +
  geom_hline(yintercept=.33,lty=2) + 
  geom_vline(xintercept=0,lty=3) + 
  facet_wrap(~item)
  scale_x_continuous(limits=c(-2,3),expand = c(0,0), name = "Time (s)") + 
  scale_y_continuous(limits=c(0,1),expand = c(0,0), name = "Proportion correct looking")

## 3. BY PARTICIPANT ANALYSIS
# This will look even worse... 
# 1. Aggregate your data by time, trial type, and item.
ms <- aggregate(correct ~ t.crit.binned + trial.type + subid, d, mean)

# 2. Plot your data as before, but with different facets for different participants.
ggplot(ms, aes(x=t.crit.binned, y=correct, colour=trial.type)) +
  geom_line() +
  geom_hline(yintercept=.33,lty=2) + 
  geom_vline(xintercept=0,lty=3) + 
  facet_wrap(~subid) +
scale_x_continuous(limits=c(-2,3),expand = c(0,0), name = "Time (s)") + 
  scale_y_continuous(limits=c(0,1),expand = c(0,0), name = "Proportion correct looking")

## 4. POSITION SALIENCE ANALYSIS
# 1. Aggregate your data by time, trial type, and item.
ms <- aggregate(correct ~ t.crit.binned + roi, d, mean)

# 2. Plot your data as before, but with different facets for different participants.
ggplot(ms, aes(x=t.crit.binned, y=correct, colour=roi)) +
  geom_line() +
  geom_hline(yintercept=.33,lty=2) + 
  geom_vline(xintercept=0,lty=3) + 
  scale_x_continuous(limits=c(-2,3),expand = c(0,0), name = "Time (s)") + 
  scale_y_continuous(limits=c(0,1),expand = c(0,0), name = "Proportion correct looking")

## 5. BLOCK ANALYSIS 
# Add a variable coding first vs second half
d$half = as.factor(ifelse(d$stimulus <= 9, 1, 2))
# 1. Aggregate your data by time, trial type, and half
ms <- aggregate(correct ~ t.crit.binned + trial.type + half, d, mean)

# 2. Plot your data as before, but with different facets for different participants.
ggplot(ms, aes(x=t.crit.binned, y=correct, colour=trial.type,linetype=half)) +
  geom_line(size=2) +
  geom_hline(yintercept=.33,lty=2) + 
  geom_vline(xintercept=0,lty=3) + 
  scale_x_continuous(limits=c(-2,3),expand = c(0,0), name = "Time (s)") + 
  scale_y_continuous(limits=c(0,1),expand = c(0,0), name = "Proportion correct looking")

## 6. DWELL TIME IN WINDOW ANALYSIS
# 1a. Compute the average dwell time on the target (correct) referent in a timw window. Which window should we choose? Aggregate correct looks by trial type and subject in that window; compute error bars.
window <- c(.2,2.5)
mss <- aggregate(correct ~ trial.type + subid, 
                 subset(d,t.crit.binned > window[1] & t.crit.binned < window[2]), 
                 mean)
ms <- aggregate(correct ~ trial.type, mss, mean)
ms$cih <- aggregate(correct ~ trial.type, mss, ci.high)$correct
ms$cil <- aggregate(correct ~ trial.type, mss, ci.low)$correct

# 1b. Make a bar graph of mean dwell times. Tip: use geom_bar()
ggplot(ms, aes(x=trial.type, y=correct, fill=trial.type)) +
  geom_bar(stat="identity", ylim = c(0,1)) +
  ylab("Proportion correct looking") + 
  geom_hline(yintercept=.33,lty=2) + 
  geom_errorbar(aes(ymin=correct-cil,ymax=correct+cih,width=.2))

# 2a. Compute the average dwell time on the target (correct) referent in pre/post word timew window. Which windows? Aggregate correct looks by trial type and subject in that window; compute error bars.
window1 <- c(-.8,.2)
window2 <- c(.2,2.5)
d$twindow = as.factor(ifelse(d$t.crit.binned <= window1[2], 1, 2))
mss <- aggregate(correct ~ trial.type + subid + twindow, 
                 subset(d,t.crit.binned > window1[1] & t.crit.binned < window2[2]), 
                 mean)

ms <- aggregate(correct ~ trial.type + twindow, mss, mean)
ms$cih <- aggregate(correct ~ trial.type + twindow, mss, ci.high)$correct
ms$cil <- aggregate(correct ~ trial.type + twindow, mss, ci.low)$correct

# 2b. Make a bar graph of mean dwell times by window.
dodge = position_dodge(.9)
ggplot(ms, aes(x=trial.type, y=correct, fill=twindow)) +
  geom_bar(stat="identity",position=dodge, ylim = c(0,1)) +
  ylab("Proportion correct looking") + 
  geom_hline(yintercept=.33,lty=2) + 
  geom_hline(yintercept=.5,lty=2) +   
  geom_errorbar(aes(ymin=correct-cil,ymax=correct+cih,width=.2),position=dodge) 

